Example Program
Global Alignments
Computing an optimal global alignment between two sequences.
File "alignment.cpp"
1#include <iostream>
2#include <seqan/align.h>
3#include <seqan/graph_align.h>
4
5using namespace std;
6using namespace seqan;
7
8int main()
9{
Here are some sequences to align:
10    DnaString seq1 = "atcgaatgcgga";
11    DnaString seq2 = "actcgttgca";
Now we choose a scoring scheme with affine gap costs ("gap open" == -2, "gap extend" == -1).
12    Score<int> score(0, -1, -1, -2);
Example 1: We use Align to align the two sequences. Since we does not specify an algorithm tag when we call globalAlignment, a suitable algorithm (Gotoh) is automatically choosen.
13    Align<DnaString, ArrayGaps> align;
14    resize(rows(align), 2);
15    assignSource(row(align, 0), seq1);
16    assignSource(row(align, 1), seq2);
17
18    cout << "Score = " << globalAlignment(align, score) << endl;
19    cout << align << endl;
Example 2: We now choose explicitely the algorithm MyersHirschberg. Since this algorithm always works on Levenshtein distance, score is ignored here. Therefore, this algorithm computes a different alignment and returns a different score.
20    cout << "Score = " << globalAlignment(align, score, MyersHirschberg()) << endl;
21    cout << align << endl;
Example 3: We now do the same as in case 1, but now we use an Alignment Graph for storing the alignment. Here we use Hirschberg's algorithm.
22    typedef StringSet<DnaString, Dependent<> > TStringSet;
23    typedef Graph<Alignment<TStringSet, void> > TAlignmentGraph;
24
25    TStringSet string_set;
26    appendValue(string_set, seq1);
27    appendValue(string_set, seq2);
28    TAlignmentGraph alignment_graph(string_set);
29
30    cout << "Score = " << globalAlignment(alignment_graph, score, Hirschberg()) << endl;
31    cout << alignment_graph << endl;
32
33    return 0;
34}
SeqAn - Sequence Analysis Library - www.seqan.de